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ABSTRACT 


The results of research relating to the feasibility of using 

* 

a low gravity environment to model geophysical flows are presented 
in this report. Atmospheric and solid earth flows are considered. 
Possible experiments and their required apparatus are suggested. 
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INTRODUCTION 

The advert of the space shuttle presents exciting possibilities for 
developing physical models of geophysical flow phenomena that cannot be 
adequately modeled in the earth's gravitational environment. 

During the past decade significant advances in the development of 
theoretical models of geophysical phenomena have taken place. For example, 
the theory of new global tectonics has led to the development of many new 
mathematical models of the solid earth. These have had a significant 
bearing on the understanding of such important phenomena as earthquakes. 

I 

In addition to the models associated with the solid earth, an atmospheric 
circulation model has also been recently developed. Many of the components 
of these models could be validated if a zero G environment were available. 

It is known that dynamical processes in various parts of the earth 
are responsible for variations in the length of the day. These variations 
comprise three distinct components: (1) seasonal fluctuations on the order 

of 1 X 10“^ sec., (2) irregular decade fluctuations on the order of 5 x 10“® 
sec., and (3) a secular increase in the length of the day by about 1 x 10“^ 
sec. per century. The secular increase is associated with angular momentum 
transfer of the earth to the moon caused by the action of gravitational 
torques associated with the tidal bulge. Seasonal fluctuations are caused 
by torques on the mantle produced by the combined effect of atmospheric 
winds and ocean currents. The amplitude of the decade fluctuations is too 
large to be accounted for in terms of interactions of the ocean and atmosphere 
and geophysicists generally agree that these fluctuations must therefore be 
due to angular momentum transfer between the mantle and the liquid core. 

The nature of the stresses that couple the core to the mantle must account 
for the fluctuating torques at the core mantle interface which are implied 
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by the decade fluctuations. The specific nature of these stresses cannot 
be determined without detailed theoretical calculations of specific models 
of the coupling process. In the past, symmetric models have typically 
been considered but there now exists substantial and growing evidence— 
to which satellite observations have made a significant contribution— 
that render symmetric models increasingly inadequate and which demand 
refinements. These refinements must reflect dynamical processes within 
the earth if they are to provide the keys to the earth's past and future 
evolution. 

In addition to the above-mentioned problem in geodesy, the following 
problem has received considerable attention. 

Fluid motion in the liquid core of the earth is widely accepted as 
the cause of the earth's magnetic field through a dynamo action. For 
twenty years no general agreement on the driving mechanisms of this fluid 
motion has been reached, however. Both processional flow and thermal con- 
vection have been proposed and challenged as possible driving mechanisms 
of the geodynamo. The concept for an experiment described in this work 
would provide a better understanding of the problem. 

The proposed apparatus would consist of a concentric inner sphere 
and a slightly elliptical outer shell which would be made to rotate and 
precess. A dielectric fluid would be trapped in the annulus between the 
two shells and a temperature gradient would be imposed across the annulus. 
An alternating electric potential between the two shells would create a 
facsimile gravity field in the annulus. The facsimile gravity is shown 

c 

to vary as 1/r . 

The feasibility of the experiment is discussed both in terms of its 
power requirements and the differences in flow produced by the strong radial 
dependence of the facsimile gravity compared to terrestrial gravity. The 


working fluid is modelled as a constant viscosity, Boussinesq fluid and 
the characteristic value problem describing the onset of thermal convec- 
tion is derived from linearized marginal stability equations. Solution 
of the characteristic value problem shows that, for an apparatus whose 
outer shell has major and minor radii of 25 cm and 24 cm, respectively, 
with a spherical core of 10 cm radius, a 10.5 KV potential is required to 
create convection at Rayleigh numbers one order of magnitude larger than 
the critical Rayleigh number of the non-rotating case. The power require- 
ment to generate the electric field is negligibly small in comparison to 
the heating power requirements, which are estimated to be 1 watt. It is 
found that the strong radial dependence of the facsimile gravity only affects 
the magnitude of the critical Rayleigh number, but does not influence the 
mode of convection. The mode is found to be identical to the mode for 
constant gravity and for gravity varying as 1/r . At a rotation rate of 
20 rpm, it is calculated from the empirical equation of Malkus (1968) that 
precession rates less than 1 rpm would create unstable precessional flow 
in the apparatus. The experiment is to be considered feasible in any zero 
or low gravity laboratory which can provide these power requirements and 
operating conditions. 

The first part of this report describes in detail the proposed zero- 
gravity geodynamo experiment. 

The final sections of the report contain presentations regarding 
models and proposed experiments for atmospheric flow phenomena. 


INTRODUCTION 


Zero-gravity laboratories such as the NASA space shuttle offer the 
unique opportunity to construct physical simulations of three-dimensional 
and large-scale planetary flows. A currently viable area of planetary 
fluid dynamics research is the study of fluid motion in rotatings spheri- 
cal annuli in response to a variety of driving forces. The results of 
this research are useful in understanding solar rotation phenomena , motion 
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in the earth's oceans and atmosphere, planetary dynamos, and more general- 
ly, fluid motion within many planetary interiors. The objective of the 
experiment described in this work is to study the response of the liquid 
core of the earth to driving forces created both by thermal buoyancy 
effects and by the precession of the earth (see Figure 1), although the 
concept of the experiment may have broader applications. To create a 
radial, facsimile gravity field in the experiment, a near zero gravity 
laboratory is required. 

The goal of the experiment is to help resolve the apparent dilemma 
created by Higgins' and Kennedy's 'core paradox' which requires that the 
radial fluid motion necessary for the geodynamo (Busse, 1975a) must 
occur in a liquid annulus that is for the most part, thermally stably 
stratified (Kennedy and Higgins, 1973). On one hand, vigorous radial 
fluid motion in the earth’s liquid core is needed to explain the existence 
of the earth's magnetic field, while on the other hand, the vigor and 
possibly the very existence of this fluid motion is limited by the 
stable stratification of the core. 

From the amount of controversy over the driving mechanism of the geo- 
dynamo (see Rochester, et al., 1975, and Busse, 1975b), it appears that 
even the hydrodynamic flow processes occurring within rotating and 



Figure 1. Precesion and Structure of the Earth. Shaded region is 
liquid core of the Earth, with inner radius R 2 of 1300 km and outer 
radius R of 3500 km. The axis of rotation of the Earth processes 
with a period of 25, 800 years. The angle between the axis of rota- 
tion and the axis of precession is 23 1/2 ° (Malkus, 1968), 




precessing spheroidal annuli are not well understood. Understanding of 
the magnetohydrodynamic flow which must actually exist in the geodynamo 
is probably being delayed by the absence of a strong physical background 
for evaluating the coupled effects of processional and thermal buoyancy 
forces in hydrodynamic flows. 

If it is eventually determined that convection can occur in the 
liquid core, there are still strong reasons to expect the influence of 
precession to be important in determining the resulting flow patterns. 

If thermal convection is found to be inadmissable as a driving mechanism 
of the geodynamo, further study of processional influences on the flow 
patterns in the liquid core will be essential (see Young, et al., 1976). 
The fact that a large fraction of the energy dissipated by the earth- 
moon system is probably accounted for by processional flow and that the 
'core paradox' appears to Impose a restriction on thermal convection in 
the core strongly supports the need for careful considerations of pro- 
cessional effects (Young, et al., 1976). 

Core convection experiments appear to be a viable topic as a zero- 
gravity experiment because (1) there seems to be a need to realistically 
examine the effect of the 'core paradox' on geodynamo models and (2) be- 
cause core convection experiments require the particular laboratory con- 
ditions presently available only in zero or low-gravity environments, the 
ability to construct radial, spherically symmetric force fields. 

Although a physical simii ',tion of a magnetohydrodynamic dynamo is 
probably impossible to construct (Jacobs, 1974) because electrical and 
fluid-dynamic processes scale differently, another need for experiments 
of general hydrodynamic flow in rotating, spheroidal annuli, arises from 
the fact that current analytical progress and numerical studies of the 
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problem have dealt only with linear models, with thin shell or other 
limiting approximations to the fully spherical, thick shell geometry re- 
quired in a study of the nonlinear flow processes occuring within the 
earth's core (see, for example, Durney, 1968a, 1968b, and 1970, and Gil- 
man, 1975). Justification of the concept of a zero-gravity geodynamo 
experiment is ample. The important question is whether or not such an 
experiment is feasible. 

In this work the feasibility of the experiment is examined by de- 
termining the conditions for instability in the fluid, both for thermal 
convection and for precessional flow. Calculation of the conditions 
needed to produce instability then yields the minimum power supply de- 
mand which the experiment imposes on the laboratory. The experiment is 
considered feasible in any laboratory which can meet this demand. 

Four coirf)i nations of driving mechanisms are probably relevant to the 
geodynamo problem: 1) flow driven by simple thermal convection in a ro- 
tating spherical annulus; 2) flow driven by precession of a slightly 
elliptic spherical annulus (see Halkus, 1968); 3) modification of ther- 
mal convection by the addition of precession; and 4) modification of 
precession driven flows by stable thermal stratification. Although the 
last possibility may be the experiment which addresses the effects of 
the 'core paradox' most directly, it is the third possibility which makes 
the largest voltage demand on the laboratory power supply. 

Therefore, the calculations presented here are made to estimate the 
conditions needed to conduct experiment 3. It is assumed that when the 
conditions necessary for the existence of themal convection exist simul- 
taneously with the conditions which create precessional instability, both 
thermal convection and precessional flow v/ill occur. 
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The effect of the strong radial dependence of the facsimile gravity 
) is also examined by comparing the wave number of the critical 
Rayleigh number for g ~ -r with that for constant gravity and for g - 
1 . 

• 

In the absence of a rigorous mathematical treatment of the geo- 
dynamo, either by analytical or numerical methods, experiments may be 
essential for interpreting geodynamo models. The hypothetical experi- 
ment described here may provide the most direct means for evaluating the 
effect of stable thermal stratification and thermal convection on the 
precessional geodynamo model. 


EXPERIMENTAL CONCEPT 


Flow in the liquid core of the Earth is characterized by the fact 
that thermal buoyancy or (in the case of Kennedy and Higgin's 'core para- 
dox') thermal restoring forces do not act in the same direction as the 
axis of rotation of the core, nor do they act only at right angles to it. 
Because the direction of thermally-induced body forces varies and because 
the solid inner core has a radius only 0.4 times the radius of the liquid * 
core surrounding it, fluid flow in the liquid region of the core can be 
expected to be strongly three-dimensional. The ease with which this po- 
tentially complicated flow can be studied by direct observation is a pri- 
mary impetus for the development of an experiment to model the hydrodyna- 
mics of the Earth's liquid core. 

The experiment must include the dominant spherical symmetry of the 
gravitational field, the effects of thermal buoyancy and of inertial 
forces acting on the fluid. A facsimile, radially symmetric gravitational 
field can be generated by an electric field acting on a dielectric liquid 
(see Hart, 1976; Chandra and Smylie, 1972; Gross, 1967; and Smylie, 1966). 

The apparatus for a zero gravity, hydrodynamic geodynamo experiment 
would consist of a concentric inner sphere of radius R. and a slightly 
elliptical outer shell of mean radius Rq. The apparatus would rotate 
and precess, as shown in Figure 2. Trapped in the annulus between the two 
shells would be a dielectric fluid (e.g., silicone oil). A temperature 

Vor practical purposes, R^ can be taken as the major axis radius of 
the outer shell, so long as the ellipticity of the shell is small. 
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Figure 2. Schematic Diagram of Proposed Apparatus. Proposed laboratory 
apparatus consists of concentric outer elliptical shell of mean radius 
Rq, an inner sphere of radius R. . Temperature of inner annulus surface 
is T^. , of outer surface T^. Silicone oil is contained in annulus and 
the whole apparatus, including ovserving instruments is made to rotate 
at a rate u and precess at a rate 
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contrast across the annulus gap could be created by circulating a heating 
or cooling fluid within the inner core and heating or cooling the outside 
shell, as required. An alternating electric potential ±V is maintained 
between the outer shell and the core to produce a facsimile gravity field 
in the model. 

A shaft (not shown in the figure) supports the core inside the shell 
and provides access for heating and cooling the inner core as well as for 
temperature measuring instrumentation, electrostatic power supply and il- 
lumination for flow visualization. 

Inertial forces acting on the fluid as a result of the precession of 
the apparatus would create turbulent motions within the annulus for rates 
of rotation and precession greater than some critical values. The effect 
of precessional forces on fluids contained in rotating, precessing and 
slightly elliptical cavities is to create a cylindrical shear layer ex- 
tending between ±30'’ latitude (see Figure 3). Fluid in the central re- 
gions of the cavity has a general retrograde (westward) drift while fluid 
outside of the shear layer shows prograde motion (Malkus, 1968). This 
cylindrical shear layer will undoubtedly be modified somewhat by the 
presence of the solid inner core. 

Thermal buoyancy forces would be created in the model by the inter- 
action of density gradients and the facsimile gravity field. When unstable 
thermal gradients would be applied, the chiefly east-west flow generated 
by precession would interact with the north-south flow caused by thermal 
convection. In the case of stable thermal stratification, the radial com- 
ponents of the precessional flow would be suppressed by thermal buoyancy 
forces . 


Figure 3. Precessional Shear Layer. The cylindrical shear layer 
observed by Malkus (1968) within a fluid contained in a rotating 
and processing elliptical cavity is shown in exaggerated form in this 
figure. With increasing rates of rotation and precession the shear 
layer becomes unstable, developing wave"! ike motions and finally 
turbulence. 



Because the electric field generating the facsimile gravity field 
cannot be made strong enough to overcome the external gravity in the 
terrestrial laboratories without creating electrical breakdown of the 
dielectric fluid (Hart, 1976), the experiment must be performed in zero 
or low gravity laboratories. 

To permit visualization of flow fields by tracer motions or dye 
streaks, the upper hemisphere of the outer shell would be constructed 
of glass or plexiglass with a thin, transparent coating of a metallic 
oxide to make it electrically conducting. The inner core could also be 
constructed of coated glass or plexiglass, permitting the use of shadow 
graph or Schlieren flow visualization techniques (see Hart, 1976). 

In contrast to the flow visualization needs for convection experi- 
ments in rotating, spherical annuli, the needs for visualizations and 
data obtained from the experiments described here include making records 
of east-west fluid motion as well as north-south fluid motion. The mo- 
tion on latitudinal planes is important, as it reflects the contribution 
of precessional flow to the resulting fluid motion, while motion on 
longtitudinal planes Indicates the contribution of thermal convection. 

Stable thermal stratification is created by heating the outer shell 
while cooling the inner core. Unstable thermal stratification is created 
by heating the inner core and cooling the outer boundary as in Chandra 
and SmylieC 1972), Precessional instability is created by increasing the 
rate of rotation or rate of precession of the apparatus. 

The electric thermal buoyancy forces are created in the experiment 
by the interaction of the temperature-dependent dielectric constant and 
the Imposed electric field. When the dielectric constant decreases with 
Increasing temperature (see Figure 4), warm liquid seeks regions of less 
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intense electric field while cold liquid seeks regions of more intense 
electric field (Chandra and Smylie, 1972). 

Another flow may be generated in the fluid due to the migration of 
free charges. This 'streaming flow' limits the accuracy of the simula- 
tion as it introduces a transport mechanism which is not found in ter- 
restrial thermal convection. An alternating electric field must be 
applied to prevent the occurrence of the 'streaming flow.* 

At this point the concept of a zero-gravity, hydrodynamic geo- 
dynamo experiment has been defined. le physical design of the appara- 
tus must wait until the needed boundary conditions have been determined, 
however. In the next section the boundary conditions for thermal con- 
vection are estimated; the boundary conditions for precessional flow 
are estimated in a later section. 


ONSET OF THERMAL CONVECTION 


The onset of convection in the annulus is determined by the criti- 
cal Rayleigh number, which is defined as the product of the Prandtl and 
Grashof numbers. The Grashof number represents the ratio of buoyant to 
viscous forces in the fluid. The Prandtl number relates temperature 
and velocity distributions in the fluid. 

An eigenvalue problem is formulated from the governing equations in 
which the characteristic value is the Rayleigh number and the mode rep- 

V 

resents the temperature distribution. The problem is formulated by 
substituting simple forms of perturbations of the state variables into 
the momentum, heat and continuity equations, yielding the perturbation 
equations. An exponential time-dependence of these perturbations is 
then assumed, and the marginal stability equations are derived from the 
perturbation equations. The eigenvalue problem is formulated directly 
from the marginal stability equations. (See Chandrasekhar, 1961.) 

In this work, derivation of the perturbation equations is taken 
directly from the work of Durney (1968a). The derivation is only 
briefly described here to provide the background needed to understand 
the marginal stability equations. 

In calculating critical Rayleigh numbers by this method, it has been 
assumed that the apparatus consists of two non-rotating, concentric 
spheres. . The first assumption is made in order to decouple the different 
modes of convection, simplifying the problem by eliminating terms in the 
governing equations (see Durney, 1968b). As a consequence, the calculated 
Rayleigh number represents a minimum Rayleigh number of interest in the 
experiment. Fluid viscosity is assumed to be constant. 
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The second assumption is reasonable since the inner boundary is 
spherical and the outer boundary is nearly spherical. 

Before deriving the eigenvalue problem, we must determine the form 
of the electrical "gravity" field. The next section is devoted to this. 

A, Electric Body Forces 

The electric bodly force per unit volume exerted on the fluid is 
(Chandra and Smylie, 1972) 

f = 1/2 P'v [e/ (H) 3 - 1/2 e/ (I f) vr (1) 

I ,0 

P.O 

where primes indicate flow-induced quantities and subscripted zeroes 

indicate stationary values. is the electric field strength, e is 

the dielectric permittivity, T is the temperature. 

The permittivity of a material is calculated by multiplying the 

material's dielectric constant and the universal constant e^, the per- 

-12 

mittivity of free space. In mks units has the value 8.854 x 10 
farads per meter. 

The dielectric constants of a material is defined as the ratio 
of the electric field strength in a vacuum to that in the material, for 
the same distribution of charge (Smyth, 1955). Another definition makes 
use of the ratio of the capacitance of a flat plate condenser with a 
vacuum between the plates (C^), and the same condenser with the material 
between the plates. The dimensionless dielectric constant is defined as 
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A typical dielectric constant for silicone oils used in convection ex- 
periments is K = 2.65 Csee, for example, Chandra and Smylie, 1972). 

The departure of density from its stationary value is 

p'=-aPj' (2) 

and the fluid state is given fay 

T = r+T^ 
p = p’ Po 

P = p' + PQ 

where p is the pressure and a is the volume coefficient of expansion. 
Again, subscripted zeroes indicate stationary conditions (no fluid 
motion, purely hydrostatic pressure field, temperature distribution 
that of pure conduction), and primes denote flow-induced quantities. 

Substituting (2) into (1) and dividing by gives the electrical 
body force per unit mass 

F ■ -1/2 ar V [e/(||.) -1/2 1/p^ e/ vr 

Assuming density and permittivity changes are small, (|f* ), ^ will be 

dp i 

independent of the spatial coordinates and the body force per unit mass 
becomes 

(t) T,o < - V2 1/P„ e„2 (II) vr 

J P jO 

It is the curl of the body force which we will use in the equation 
of motion, or 
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VxF = -l/2o( — ) VT'xVE^ -V2 1 /p„(|t) Ve/ x VT' 

''3p''T,o ° ® 

Equivalently, since the two vector components of V x F are co-1 inear 


7 xt= 1/2 tl/,P„ (-|f) 


p,0 




The permittivity coefficients must also satisfy the relation 



from thermodynamics (Chandra and Smylie, 1972). Thus the curl of the 
body force may be written 


7 (If) ’T' X 7 e/ 

° p. 


(3) 


By analogy to the terrestrial thermal buoyancy force 
F =«T’ V4> 

where g =v<|) and <|) is the geopotential, it can be shown that the elec- 
tric, facsimile gravity in (3) may be written 

2oPo 

The electric field E in a spherical capacitor with inner radius R-, 
outer radius R^, and annulus gap a is (Moore, 1973) 
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where V is the voltage across the capacitor. Inserting this in (4) we 
obtain the useful expression for the electric facsimile gravity 



The radial dependence of remains a major difference between the 
facsimile gravity field and the gravity field which actually exists with- 
in the interior of the earth, although this difference may be most no- 
ticeable only at the onset of thermal convection (Gilman, 1976). Other 
kinds of radial dependence can be produced in facsimile gravity fields 
by the use of other geometries for generating the electric field, as 
shown in Table I, 



Table I 


Geometry 

Electric Field 

Facsimile Gravity 

plane 

constant 

none 

cylindrical 

1/r 

l/r^ 

spherical 

1/r^- 

1/r= 


Only with spherical geometry does the facsimile gravity field posess 
the spherical symmetry found in the earth's gravitational field, however. 

B. Formulation of the Eigenvalue Problem 

Consider a stationary spherical annulus of thickness a filled with 
a fluid of density p^. The acceleration of gravity within the annulus 
is 
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The outer radius of the annulus is R^, its inner radius tiRq = The 
velocity field u, temperature T, radial coordinate r and time t are 
scaled by the following definitions: 

u = ” u‘j T « |Tq| T'; r = R^r' 

R 2 

t = — t'; gCr) = g(R ) g'(r‘) 

K ° 

K is the thermal diffusivity and T^ is the negative temperature of the 
outer boundary of the annulus. (The temperature at the inner boundary 
is assumed to be zero.) All primed quantities are dimensionless. 

Dropping the primes in (6), the nondimensional momentum, continuity 
and heat equations of the problem can be written as (Durney, 1968a)^ 

— — 7 X u - V X 7^ u = -T V X (u *V) u + R,V X g(r) r T 

at > ' 

V ' u = 0 
and 

t It - V^) T = - V(u T) 

Pj^ is the Prandtl number 

%t should be noted that this is the curl of the momentum equation, as 
the pressure term Vp is absent. 
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P.. = 


K 


and is the Rayleigh number based on the radius of the outer annulus 
boundary, 

3 


R, = ^iTn I g(Rn) Rn 
K V 


Another common form of the Rayleigh number is based on the annulus 
depth a , 

oITgl 9(Ro> 




K V 


and is related to the outer radius Rayleigh number R-j by (Durney, 1968a) 

R, = R, C 4- ] ^ 


Equivalently, 


Ra = Ri (1- n)' 


Substituting equation (5) into the definition of R^ we obtain 


ATV^( If- ) 


Hl 


R-" Ro (i-n)2 




K V 


the useful form of the Rayleigh number. 

Durney (1968a) has derived the perturbation equations fay assuming 
a temperature distribution of the form 

T ■ — [ ?■ - 1] iKr.t) -I- 6(r, t) 

1-t) 


I 


The first term represents the temperature field of pure conduction with 

boundary conditions T = 0 at the Inner boundary and T = -1 at the outer 

2 

boundary. It Is a solution of Laplace's equation V T = 0 In spherical 
co-ordinates. t|» (r,t) represents the mean distortion of the temperature 
field by convection and Is a function of the radial coordinate 

A 

and time, e (r, t) Is expanded In spherical harmonics by writing 

0 (r,t) = Z 0L ^ (ir,t) (e.tjj) 

L»ni ’ 

where yj^”* is the spherical harmonic and Is a known function of e , , 

and wave numbers L and m. Thus the temperature field Is completely de- 
termined by specifying only two functions, 0 , _ C^jt) and i|»(r,t), both 
of which are functions of radius and time only. 

Similarly a specific form for the velocity is assumed using the 
pololdal vectorial _ 

U 9III 

which has the following components in spherical co-ordinates 


P (»') = 


t~- Pl.,m t)yL"’(0,<}.) 


ap, „ (r,t) y, (e,<p) 


I ZUR 

r 


1 a Pi.m 


r sin 0 
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Thus the velocity field is completely described by determining a function 

Pj_ ^ (r,t) which is a function of radius and time only* 

The convection problem is solved by finding the required functions 

P, (r*t), 0, „ (r,t) and tj/(r,t). From the momentum, continuity, and 
L,m 

heat equations, Durney (1968a) derived three equations governing these 
variables, the perturbation equations 


D, 


Pi _ “ 

L,m 1 


g(r) 0 


L,m 


8(a) 






1 

4ifr^ 


Z (L+1)L 
L,m 


ar 



a 


L,m 


) 


8(b) 


and 


36 


L,m 


at 


- D, 0 


L,m 


(L-H)L 

r 


P. m C 
L,m *• 


TT^ 


_ M 1 
ar *' 


where 0|_ is the differential operator 



dr 


iktilk 


8(c) 


and, for his own convenience, Durney has redefined P, as 

L gill 


P 


L,m 


P 


L,m 


(r,t) 


r 
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The boundary conditions on the temperature variables are tp = 

= 0 at r = n> !• Either rigid or free boundary conditions are used for 
the velocity, as required* 

Since we are interested in those small perturbations which start 
convection, 0, Ui » and p, „ may be assumed to be small. By ignoring 
products of the perturbations with themselves and their derivatives the 
linearized perturbation equations may be obtained 


2 

“l \.n, 

= R-j g(r) 

®L,m 


9(a) 

It ■ “l'I' 

(L+1)L 

■ 

i|> = 0 


9(b) 

^®L,m 

■ »L »L.n, = 

(L+DL 

P ^ 


3t 


(1-n) 

9(c) 


Note that different values of wave number L are decoupled and that 
the equations are m independent. Since the term containing a sum 
over m terms in 8{b) has vanished in the linearization process, how- 
ever, the governing equations themselves are independent of m and the 
subscript will be dropped. 

The marginal stability equations are obtained from the linearized 
perturbation equations by assuming an exponential time dependence in 
the state variables. The perturbed temperature distribution will be 
described, for example, by the product of a function of time and a func 
tion A describing the amplitude of the perturbation in terms of the 
spatial coordinates: 


The conditions for stability are clearly 


p > 0: unstable 
p = 0: marginally stable 
p < 0; stable 

Thus the marginal stability equations are obtained by setting the ex- 
ponent p to zero in the perturbation equations. 

The three functions used to describe the temperature and velocity 
fields are rewritten with the exponential time dependence 

0L (r,t) = eP^ e^Cr) 

i(j(r,t) - eP^ ifi(r) 
and 

P(_(r,t) = eP* p^Cr) 

Substituting these definitions into the linearized perturbation equations 
under marginal stability cond-itlons {p=0)j all time derivatives dis- 
appear and the time-dependent teni;s become unity. The result is the 

10(a) 

10(b) 


marginal stability aquatic s 

Pl *^1 

D 1^ + 4, = 0 

r 

D, 0, = P, — 0— 

^ r^ ^ (1-n) 


10(c) 
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Equation (10c) can be written as 


P = idiB: 


»C ® L 


where (1+1)1.. Inserting this expression into equation (10a) yields 
the eigenvalue problem in terms of the temperature variable e^, 

D|_ = R, g(r) e L 


in which the Rayleigh number Ri is the characteristic value and the 

j * 

temperature distribution is the corresponding mode. By operating 
the differential operator Dj^ upon itself the differential operator 
may be obtained 


2 _ d" 


>• dr^ 


?h_ d' 

d? 


. h 


; ^ note that the definition of the non-dimensional gravity term 

--1 '-j 

(r) when applied to the electrical gravity yields (see equation (6)) 


g'(r) 


9%) r' 


Using the definitions of Dj_ and the eigenvalue problem may be 


written as 


, 18 dSe, 

— _ — , , — 


96-3L, 


(continued) 






dkiginac page is 

?OOR QUALITY 
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168 - 24L 


1 




0L^ 

dr^ 


3L^ -42L^+72 




6L^ •* 12L^ 


d0 

L 

dr 


- 2l2 

H ^'*1 


0. = 


nL-| g 

O"' L 

n-i r 


( 11 ) 


This is the basic equation of the eigenvalue problem. With the six 
boundary conditions derived in the next section, foniiulation of the 
problem describing the onset of convection is complete. 

No-slip (rigid) boundary conditions on the velocity and constant 
temperature at inner and outer annul ui boundaries are used. The tempera- 
ture boundary condition is 


0|^ (r) = 0 at r - n. 1 

and 

if((r) =0 at r = n, 1 


(12a) 


The no-slip condition requires that all components of velocity 
vanish at the boundaries, or 

p^(ir) ^ p^( 0 ) = P|^ (<)>) = 0 at r = n. 1 . 

By the definitions of the components of the velocity, it is apparent that 
only two conditions are needed to make u vanish at the boundaries. These 


are 
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Pl^(r) = 0 at r = n, 1 


and 


ap^Cr) 

dr 


= 0 at r = n» 1 


Using the former in equation (9c) we obtain D. 6. = 0 at r 1. (12b) 

P. 

Substituting the definition P^ = — into (10c) 

L, 


JL 


"" y,^ i 

r n-j 


Pl 


The requirement 


dPL(*') 

dr 


0 can be restated using this expression as 


W Cr“ 6|_] = 0 


(12c) 


Upon substituting the definitions of the differential operators in- 
to equations (12) we obtain the boundary conditions in their useful 
form: 


0^ = 0 at r= n> 1 


(13a) 


d^ei ^ 2 d8|_ 


dr" 


+ - 

»’ dr 


9 = 0 at r=njl (13b) 


and 


d^e, 


dr" 


d^ 0j 

dr^ 


6-U d 0, 


dr 


(13c) 


2U 


0, = 0 at r = n, 1 


Thus equation (11) is the governing equation for the eigenvalue 
problem and equations (13) are the six boundary conditions needed to 
solve it. 

C. Numerical Solution 

The eigenvalue problem (11) and associated boundary conditions (13) 
were solved by using the finite-difference method and the edition 5 Inter- 
national Mathematics and Statistics Library (IMSL) eigenvalue subroutine 
EIGZF (see Appendix I). The domain r s 1 was modelled with 37 nodal 
points. It should be noted that the only solution to (10b) which 
satisfies the boundary conditions is i/j (r) = 0. The temperature distri- 
bution within the annulus is therefore given by 

T = —[1 -1] + 6(?, t) 

1- n 

Because (11) is homogeneous, any multiple of a given solution 0L(r) is 
also a solution. 

Solutions were obtained for n = 0.4, rigid boundaries and g - 

r^ 

to estimate critical Rayleigh numbers for the experimental apparatus. 

To evaluate the influence of the strong radial dependence of the radial 
gravity on the mode of convection, solutions were obtained v/ith free 

C 

boundaries (see Appendix III) with n = 0.8, g - 1/r . These results were 

compared with the results of Durney (1968a) and Gilman (1975) which 

2 

were obtained with constant gravity and with g - 1/r , respectively. 
Verification of equation (11) as a model of thermal convection in the 
annulus was obtained by (1) repeating the result of both Durney (1968a) 
and Gilman (1975) and (2) demonstrating that the model predicts higher 


- 


I 
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critical Rayleigh numbers for rigid boundaries than for free boundaries. 
The results are summarized in Table II. 


SoTn. 

Gravity 

a 

Boundaries 

L 

5a 

1 

const. 

o 

• 

CO 

free 

9 

740 

2 

1/r^ 

0,8 

free 

9 

596 

3 

1/r^ 

0.8 

free 

9 

426 

4 

1/r^ 

0.4 

free 

3 

142 

5 

1/r^ 

0.4 

rigid 

3 

295 


Table II 

Summary of Numerical Results 

Solutions (1) and (2) reproduce within 5 % the results of Durney 
(1968a) and Gilman (1975), respectively, indicating that the linear 
model of convection represented by equation (11) with the free surface 
boundary conditions is in general agreement with nonlinear models of 
convection. Solutions (4) and (5) demonstrate that the linear model 
with rigid boundary conditions is also consistens with physical Intui- 
tion (Durney, 1976), because the model predicts higher critical Rayleigh 
numbers for rigid boundaries than for free boundaries. Higher Rayleigh 
numbers for rigid boundaries are to be expected, since in that case 
thermal buoyancy forces must overcome viscous forces both in the bulk 
of the fluid and at the fluid boundaries for convection to occur. With 
free boundaries buoyancy forces must overcome viscous forces only in 
the bulk of the fluid to initiate convection. 


Comparison of solutions (1), (2), and (3) shows that the effect of 
the strong radial dependence of the facsimile gravity only alters the 

5 

value of the critical Rayleigh number. For g 1/r , the critical Ray- 
leigh number occurs at a wave number L = 9 (see Figure 5), as did the 
critical Rayleigh numbers calculated by both Durney (1968a) and Gilman 
(1975) for constant gravity and g 1/r , respectively. The spectrum of 
Rayleigh numbers over a wide range of wave numbers is shovm in Figure 5. 
Figure 6 shows the critical mode and temperature distribution in the 
rigid boundary case. 

Since the effect of rotation on thermal convection in spherical 
annuli is to suppress fluid motion and therefore increase the critical 
Rayleigh number (Gilman, 1975), the results calculated here represent 
the bottom of the range of important Rayleigh numbers for the experi- 
ment. It will probably be desirable to conduct experiments over the 
range of Rayleigh numbers from the minimum to at least an order of 
magnitude or more above the minimum. 




boundaries 



Wave Number L 


Figure 5a. Rayleigh Number vs. Wave Number L for n = 0*4 and rigid 
boundary conditions. Minimum Rayleigh number of 295 occurs at a 


wave number of 3. 
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i 


Wave Number L 


Figure 5b. Rayleigh number vs, wave number L for n = 0.8 and free 
surface boundary conditions. Minimum Rayleigh number of 426 occurs 
at a wave number of L = g and L = 10. 
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fi;iD ivs 


Figure 6a, Critical Mode. Plot shows e L » 3, from inner annulus 
boundary to outer annulus boundary, 6j_ is arbitrarily normalized so 
that its largest element is unity. For a free surface boundary solu- 
tion to the same problem, see Appendix III. 


l.ttt 







Precessional Flow 


The equation of motion of a viscous fluid inside the spheroidal 
cavity of a precessing rigid body which is rotating about its axis with 
a constant angular velocity w^ is 

2 S X q + q • 7q = -7P + E 5 

where is the angular velocity of the rotating frame of reference, p is 
the pressure and q is the fluid velocity {Busse, 1968). 

This equation illustrates the fact that equality of the Eckman num- 
bers 



for two precessional flows will establish dynamic similarity between, for 
example, the liquid core of the earth and a laboratory experiment. Some 
physical properties of and vertical velocity in the earth's liquid core 
are shown in Table III. These figures yield an Eckman number for the 
core of 

E = 6.72 X 10“^^ 

based on the radius of the earth's outer core, 3500 km. Figure 7 shows 
that in a feasibly-sized apparatus rotated at speeds below 5000 rpm, Eck- 
man numbers in a laboratory experiment will be much larger than in the 
earth's core. Although it appears impossible to establish dynamic 
similarity between the experiment and the core, previous experiments with 
fluid motion in rotating and precessing cavities indicate that instabili- 
ties and turbulence can be made to occur for E >10"^ (Malkus, 1968). Thus 
fluid motion of interest can be created in a physical simulation of the 
earth's liquid core. 




• R (meters) 


Figure 7. Eckman numbers in the laboratory experiment. R is the mean 
radius of the outer shell, (v = 12 CS, after Chandra and Smylie, 1972} 
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Property 

Value 

Reference 

Specific Heat 

7.12 X 102 Jkg-i®K-i 

Frazer (1973) 

Coefficient of Thermal 
Expansion 

4.5 X 10-^®K“i 

Frazer (1973) 

Thermal Conductivity 

60 W 

Frazer (1973) 

Radial Fluid Velocity 

3 X 10"'* m s-i 

Frazer (1973) 

Kinematic Viscosity 

0.6 cent i stokes 

Gans (1972) 

Mass Density 

13 gm cm“3 

Gans (1972) 


Table 3. Physical Properties of and vertical 
velocity in the Earth's liquid core 
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Onset of instability In a fluid contained ih a rotating and pro- 
cessing spheroid can be determined from an empirical equation relating 
the ratio A of the maximum toroidal velocity in the fluid to the speed 
of the periphery of the container. 

Instability occurs when 

A > (5.0 ± 0.3) E (14) 

where E is the Eckman number. 


is the mean radius of the spheroid, and w is the rotation rate about 
the sphere's minor axis (Malkus, 1968). A is determined from the rota- 
tional velocity of the container (w) and the rate of precession (Jl) 

A = f f(E) ' (15) 
where e is the ellipticity of the spheroid. 


(see Malkus, 1968). Ip is the moment of inertia of the spheroid about 
Its "polar" (minor) axis and is its moment of inertia about t.-o 
"equatorial" (major axis). 

In (15) f(E) is defined by 
f(E) = 0.4, E > 10“^ 
f(E) : E < 10“^ 

Combining the definition of the Eckman number and expression (15) 


for A with equation (14), we find the required rate of precession for 
instability to be 
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For an elliptical body with a major axis R„^. and * ;nor axis R„. 

maj mm 

the moment of inertia about the pole (minor axis)-fs 

Ip “ f 


and the moment of inertia about the equator (major axis) is 

[ = ^ 

5 


^ ^'^maj^ '^min^^ 


so the dynamic el'iipticity of the body is 

2 ? 

R ^ - R ^ 

„ maj mi n 


2R, 


maj 


Malkus (1968) found that instability in the processional flow first 
manifests itself as wave motions with the cylindrical shear layer de- 
picted in Figure 3. The waves move retrograde (west) relative to the 
rotation of the spheroid. For a given rate of precession, instabilities 
are intensified by decreasing the Eckman number of the flow. 


EXPERIMENTAL APPARATUS 


The critical Rayleigh number for convection in the experimental ap- 
paratus was calculated to be R^ = 295. Using the definition of the Ray- 
leigh number we find 

2 = 295 

for convection to occur when the apparatus is not rotating. 

The dimensions of the hypothetical experiment considered here are 

listed in Table IV. 


Table IV 

Dimensions of Hypothetical Experiment 
Symbol Term Numerical Value 


n 

aspect ratio 

0.40 

'>1 

inner radius 

10 cm 

''o 

outer radius 

25 cm 

a 

annulus depth 

15 cm 

AT 

temperature difference 

5“C 

e 

ellipticity 

0.05 


The properties of the hypothetical wording fluid used in this ex- 


periment are listed in Table V. 




40 


Table V 

Properties of Hypothetical Working Fluid 


Symbol 

Term 

Numerical Value 

Reference 

a 

coef. of expan. 

1.08 X 10"^°C 

1 

(^) 

coef. of permit. 

3.32 X 10“^^ f/m-«C 

1 

K P"" 

th. diffusivity 

6.40 X 10“V/s 

2 

V 

kin. viscosity 

10“V/s 

2 

p 

density 

873 kg/ft? 

2 

"o 

permittivity const. 8.854 x lo"^^ f/m 

3 

c 

dielectric const. 

2.60 

1 

Cp 

sp. heat 

1.9 X 10^ j/kg'’C 

2 


The references are : (1) Chandra and Smylie (1972); (2) Hart (1976); (3) 
Moore (1973). 

From the definition of the Rayleigh number we find 

V = 3.34 KV 

for convection to occur in the hypothetical apparatus (n = 0.4) without 
rotation and with a 5"C temperature contrast across the annulus. We note 
that with a thin annulus gap (n = 0.8) the voltage requirement is 

V = 0.55 KV 

under the same operating conditions. To achieve Rayleigh numbers one 
order of magnitude larger than the critical Rayleigh number in the sta- 
tionary apparatus requires 


V =10.5 KV for n = 0.4 
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V "1.76 KV for n = 0.8 

The capacitance of the apparatus is given by 

C = ^ ^ gRj Rf) 
a 

which is the capacitance of a spherical condenser (Moore, 1973). For 
the hypothetical experiment 

C = 48.2 X IC"^^ farads 

The electrical energy stored in the experimental apparatus is given by 
(Moore, 1973) 

w = 1/2 CV^ 
or 

w = .003 watt-sec 
for a voltage of 10.5 KV. 

For a rotation rate of 20 rpm the Eckman number is 7.64 x 10 . 

The rate of precession required to produce instability is calculated 
from equation (16) for an inclination of 30® between rotation and pre- 
cession axes, 

= .56 rpm 

Modeling the circuit needed to create the electric field as a one 

-12 

ohm resistance and capacitance of 48 x 10 f in series, the total im' 
pedance of the circuit is 5.50 x 10^ SI . The real power supplied by a 
60 hz voltage supply of peak output is then 

P = 1 .60 10"^ watts. 

Heat flux by conduction is (Krieth, 1973) 


where k is the thermal conductivity. The fluid properties in Table II 
give a conductivity of k = k CpP = 1.06 x 10"^ 
heat flux Is thus 
q = 1.11 watts 

for a temperature contrast of 5*C. 

The total power demand on the laboratory is therefore on the order 
of 1 watt. To achieve Rayleigh numbers up to one order of magnitude 
greater than critical requires up to about 10 KV. At a rotation rate 
of 20 rpm, the empirical equation of Ma'^kus (1968) predicts that a pre- 
cession rate of 0.56 rpm is required to initiate instabilities. 


CONCLUSION 


A. Current Work 


A well-defined experiment could be essential for providing the 
physical basis to evaluate and compare geodynamo driving mechanisms. Of 
particular interest is an experiment to evaluate the effect of stable 
thermal stratification and thermal convection on the precessional geo- 
dynamo model. In this work it is found that such an experiment is 
feasible in an apparatus consisting of a concentric inner sphere of 10 
cm radius and an elliptical outer shell with major and minor radii of 
25 cm and 24 cm, respectively, provided 

(1) a potential of the order of 10 KV is maintained between 
the shells; 

(2) a temperature difference of 5°C is imposed across the 
annul us; 

and 

(3) the apparatus can be made to rotate and process at rates of 
20 rpm and 1 rpm, respectively, with a 30° angle betv/een the 
axis of rotation and the axis of precession. 

B. Recommendations fcr Future Research 

Besides the design, construction, and testing of the experimental 
apparatus important tasks which should be completed as soon as possible 
include: 

(1) Inclusion of variable viscosity and the effect of rotation in 
the numerif.al model of convection presented here; 

(2) Development of a numerical model of the precessional flow to 
verify the empirical equation of Malkus (1968) and to de- 
termine the effect of the inner core on the flow; 

(3) Combination of the two numerical models to predict the results 
- of the experiment and to verify the assumption that both con- 
vective and precessional flow instabilities will exist when 
the conditions generating them are present simultaneously. 
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Computer Algorithm for the Eigenvalue Problem 
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Eigenvalue Algorithm 

The algorithm described here solves eigenvalue problems of the form 


d ® ^ dx^ 2 3 ^^3 ^ 5 6 7 


( 1 ) 


on the domain nix <_Jl, with the six boundary conditions 


fll fl2 fl3 ^ * ^14 ^ ^ ^15^ = “ 
1' dx^ dx*^ dx*^ dx 


(2a) 


'23 


dx^ 


+ f 


24 


dx 


+ fg5y = 0 


(2b) 


and 


y - 0 


(2c) 


at X = n, A. X is the eigenvalue and coefficients c^ and f^^ are functions 
of the independent variable x and wave number w. 


The finite difference operators used to create the matrix equation 
equivalent of (1) and boundary conditions (2) are shown in Figure 8. The 
domain is modeled as a set of {n+6) nodes with n nodes on the interior of 
the domain, 4 nodes exterior to the domain and 1 node at each boundary 
(see figure belov/). 



Z 


" X ' 7* ■ " 'X jf X X 

n+1 n+2 n+3 n+4 n+5 n+6 
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At node i the finite difference equation representing (1) is 





i-2 



where h is the spacing between nodal points. Applying the finite differ- 
ence operators across the domain and applying the boundary conditions 
puts the eigenvalue problem in standard matrix form, 

[A] {y} = X [B] {y} 

in which [A] and [B] are n x n’ matrices and {y} is the n x 1 eigenvector. 
Off-diagonal elements of [B] are zero and [A] is a banded matrix with a 
bandwidth of 7. 

By using the finite erence operators the boundary conditions 
must first be stated in the useful form 

= BC (1,1) y^ + BC(1,2) yg (4a) 

and 

yj ■= BC(2,1) t BC(2,2) y^ 


(4b) 


at the left end, and 

« BC (3,1) y^ + BC (3,2) y. (4c) 

and 

y„ = BC (4,1) y^ + BC (4,2) y^ (4d) 

at the right end (see figure below). 


n 

I ^ X 

X X X H — ^ 

1 2 3 4 5 

For example, writing the second-order boundary condition in finite differ- 
ence form for the left end, we have 



6 


’i'^r 


j k 1 


m 


where the coefficients f- • are evaluaiod at x = n- 

I J 


-I 


*^^24 ^^23 


hfg^ - 2fg2 


and we have thus found 


BC (1,1) = 


hfg4 + 2f23 


^^^24 " ^^22 


Solving for y^ 


(5a) 


and 


bC (1,2) = 0.0 


(5b) 
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Writing the fourth-order boundary condition in finite -difference 
form at the left end, we have 



where again the coefficients f.. are evaluated at x = n* After solving 


for y^ and noting that 
yg = BC (2,1) y-| we have 


y^ = |(BC(2,1) + 1) 


hfig + 2f^^ 


^13 - g^ll 
*^^12 “ ^^11 


hfi2 - 2f^^ 


^5 


+ (BC(2,1) - 1) 


" ^^11 


We have now found 


2h‘f,3 - 8f„ 


2hf,2- hfj^ 


BC (1.1) = (BC (2.1) -M) jhf,2-2f,, T (^0(2.1 )-l)j { (5c) 


and 


BC (1,2) = 


hfig ■>• 2f^^ 
” ^^11 


(5d) 


Similarly it can be shown that 

BC (3,1) = 0.0 (5e) 


O 
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BC (3,2) = 


hf2^"2f23 

•^■*'24 ^’*'23 


(5f) 


BC (4,1) = 


*^^12 - ^'*'11 
^■*"12 ■*' ^‘*'11 


(5g) 


and 


BC (4,2) =(BC (3,2) + 1) ) 


«^11 - ^^13 

2f,l + hf,2 


+ (BC (3.2) -1) 


2hfl2 - h^fu 

hf,2 + 2^1 1 


(5h) 


where coefficients f^ . are evaluated at x = ji. 

Now the matrix generated by equation (3) is an n by (n + 6) matrix 
which must be transformed to an n x n matrix through the application of 
the boundary conditions. With two interior nodal points, for example, 
eight nodal points are used to model the domain: 


n 1 



1 2 3 4 5 6 7 8 


The resulting matrix equation formed by repeated application of (3) is 


1 




- 


r 1 

^11 ^12 ®13 * • * ®18 


^1 


“’ll *>12 


^4 




= X 




®21 ®22 ®23 ' • * ®28 


^2 


‘’21 '*22 j 

1 

^5) 



1 • 1 






V. J 


V 
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The boundary conditions are stated in the form 
- BC (1, 1) + BC (1, 2) yg 

yg = BC (2, 1) y^ + BC (2. 2) y^ 

y^ BC (3. 1) yg + BC (3, 2) y^ 

Yg BC (4. 1) yg + BC (4, 2) yg 

Substituting these into the matrix equation the eigenvalue problem 
becomes 



where 

a-i"^ ” a-j^ aij-j BC Osl) ®12 (B»l) ®*]‘^ BC C3»0 ^18 

®15 * ®15 **■ ®11 ®12 ■*■ ®17 ‘^18 

+ Og-j BC (1 »1 ) + agg BC (2,1) + ag^ BC (3,1) + agg BC (4,1) 

and 

a^g = agg + Og^ BC (1,2) + agg BC (2,2) + ag^ BC (3,2) + agg BC (4,2) 

At this stage the eigenvalue problem has been reduced to the form 
[A] {y} = X [B] ' {y} ’ 

which several readily available computer routines can solve. 

For convenieuce, the International Mathematics and Statistics Library 
(IMSL) edition 5 routine EIGZF was used to solve the problem, once it had 
been put in the standard form. It should be noted that the matrix equa- 
tion equivalent of equation (1) with boundary conditions (2) is not 
necessarily symmetric, requiring a sophisticated eigenvalue routine. 
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Subroistine EIGZF was found to work well for all cases and provided 
complex eigenvalues and eigenvectors. With only a few exceptions, the 
Imaginary parts of both the eigenvalues and eigenvectors were zero. 

Use of Program CONVECT 

To solve a particular eigenvalue problem the user of program CONVECT 

must 

1 . punch cards with the coefficients of the governing equation c^ 
1=1,7, and insert them In subroutine COEFF; 

2. punch cards with the coefficients f.-^-, i = 1 , 2, j = 1 , 4, and 
insert them in subroutine BCOEFF; and 

3. punch data card(s) with input values for 

ETA : left end of domain 

XF ; right end of domain 
L : wave number 

NUFiDIV: number of divisions of the domain 

The data card(s) are punched in the format shown below: 

L NUMDIV ETA XF 
(no) (110) (F10.4) (F10.4) 


The coefficients c^ are specified by assignment statements of the form 

c(l) = c^/2. 0/0X5 
c(2) C 2 /DX 4 

c(3) C 2 / 2 .0/0X3 
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and 


c(4) = 
c(5) = 
c( 6 ) = 

cC7) = 


C4/DX2 

C 5 / 2 .O/DX 


The coefficients f.. are specified by statement functions of the form 
F]](XjW), FtpCxjW), etc. 


12 ' 


where w is a wave number. Both arguments must be included, even if both 
are not used. On the printed output the program will provide, anrong other 
things, the eigenvalues and eigenvectors as complex numbers. A single 
eigenvalue or element of an eigenvector will appear as 

(real part, imaginary part) 

For the 35 x 35 matrices used to solve the critical Rayliegh number prob- 
lem, central processor times of slightly less than 20 seconds were typical. 


Verification of CONVECT 

To verify the accuracy of- the program, a sixth order eigenvalue prob- 
lem for which one solution could be determined analytically was solved. 

The verification problem was 

4+ k 4 + ^ di , 1^4 4 

dx° dx^ dx^ dx** dx'^ 

with boundary conditions 

2 3 

y = = 0 at X = 0,TT and ^ = 0 at x = 0,ir 

dx dx"^ dx 


1 


This verification problem was suggested by Dr. Bernard Durney. 
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A particular solution of the equation is 
y = Asinkx + Bcoskx 

After applying the boundary conditions the solution is reduced to 
y = Asinkx, k = 1, 2, 3, . . . n 

where k is the characteristic value and sin kx is the mode. The problem 
was put on the computer in the form 



with the boundary conditions 

y = ^ =0 at X = 0, TT 
dx*^ 

and 

4 + ^ = 0 at X = 0, u 

d^3 dx 

Thus the critical (minimum) eigenvalue was expected to be 
k = = 1.0 

and the mode y = sinx. The results are plotted in Figure 9 and show good 
convergence. The convergence of the convection problem (see Chapter III) 
is shown in Figure 10. Some instability was present when fewer than four 
interior nodal points were used to model the domain, but these instabili- 
ties disappeared as the number of equations was increased and the solution 
converged to two significant figures at 35 equations. 


value 



1«.t 




Figure 9b: Verification of CONVECT Convergence of Mode. Solid line 
is a plot of sin x for 0 < x < it and is the exact solution to the eigen 
value problem. The squares are the elements of the eigenvector com- 
puted by program CONVECT The figure shows good agreement between 
exact and computer solutions. 





s 

Number of Equations 


Figure 10. Convergence of Rayleigh Number. Solution of annulus problem 
is plotted for up to 40 equations. Some instability is present with less 
than 4 equations. 35 equations were used in the numerical solution de- 
scribed in Chapter III. 
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Program Listing and Sample Output 


61 


H 


Included in this appendix is a listing of program CONVECT and 
sample output. The output is for the rigid boundary annulus problem 
with L =!♦ n = 0.4, and NUMEQ =10. The first two pages of output 
describe the formulation of the problem, ending on the second page 
with the matrices A and B ready for input to the eigenvalue routine. 

On the third page the solution of the eigenvalue problem is printed, 
with eigenvalues representing the Rayleigh number R-| and eigenvector 
representing the ten^)erature variable e. 6 is arbitrarily normalized 
so its largest element is unity. The minimimi eigenvalue is R-| » 2063.01 
or, equivalently, R^ = (.6) » 446. The corresponding mode is 

.656 
.875 
1 .00 
.988 
.890 
.746 
.583 
.420 

If* 

.265 


.125 
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»llNfOTiCurpuT.TAKS«lKPUTtTAPE««0ur^UTt 


using the fiNtTC*ozFreft£ricE hethco* this moghah souyes > o;ie- 

DZHEN&ICNAL* SlXTH-OfiCE^ EIGENVALUE PR9iL£H 4ltH T»C EE» SitOE.%» 
T^O S£C;N3 AHO TUO rOUitTH GnOE^ BOUNOA^r COrtOITIvHS* 

PU^THEE OGCUHENTATXCN IS T3 HE POUND IN 

ccx« a.L* c 19761 PEASIBILXTV OP A Z£EO CRAVttV GE3CT6A.10 
EXPEEXPENT* HSC THESISt OEPAA7HEHT OP PEGH^NICAL 
Ei>GINEERIN6t CCLCEA30 STATE UNIVEESXTT. PSET CCLLX.NS 
OCLOEAOO 

.NOTE*- the U$C<t MUSI SUPPlT ThE SOEPPXCIENTS CP T«iE EAUATXOH 
ANO T«E ^OUNOAEV CCNOXTXONS XN SUiECUTXiiES CJEPP ANO jCGEP* t 
AESPECTIVSLT* X.NPCRTANT '/AEXA8LE3 XN ThE PEOGE-M AEE 

CN INPUT 

L HAVE NU**0EE 

NU»OXV HUHSER op OI'/XSIOnS op 50HAIN 

ETA CC8E0INATE OP lXPT ENS CP C3MAXN 

iP CCOROIKATE OP EIGHT £NO OF OCMAIN 

ON OUTPUT (TO SUSECUTINE EXGE'O 

NJHE3 NUM3EE CP EOJATXONS 

-STOtXfJt MATEXX a 3F STAftOAEO-PaP.-t EIOEnvALJE PE03L£H 
SiI,Jl ASSOCIATED JiATEIX 3 OF EISt.NVALJt P»;3Lfc'* 

curing EXECUTICH 

X POSITICK OP NODE 

-IN SPACING 3ETMEEN NODES 

CCt CCEPPICXE'NTS OP GOUEENINC E-IUATIC.n 

ECdtJ) ECUUOArT CONOITION HAT31X 

RrtSCJ RIGHT— *4N0 ilOE OF EiOATtO'l I 

ItE EEECR PARAHETEE F»Of1 EOUTI.-*E EICZP 

1 ;i*JI NUtlEC 8N INUN£0*6» HATEIA OP GC^E-t.NlNC tOUATlOM 

■(LI FEAk EEPRE3ENTATI3N OP HAVE NUH3EE ..1 


CSHNGf./::?F/ 

»ETA,xF,««CX,EL1 

ccm«c::/s:g£nv/ 

»ASTA(I5.33»(a(35»35l 


CiNENSICf. 

«A(Ll.<.ll,R.H3<2SlfaC(4t21 tCiri 


: READ AHO CHECK INPUT CATA 

* - 

35 CONTINUE 

EEA0(St2) LtMUHCIVtETAvXF 
IP tECF( 51» 1200* XO 

<tC contii.ue 


ORIGINAIJ PAGE IS 
OE POOR QGALTrXi 


WRXTf 

•JAITEI6»2) Lf KUnClVfETAtXF 


rt«irXiLIZ£ 

ioace>i 

LI\£«5 

<«£TA 

CX«(XF»ETA)/FL3ATir«UK0lVI 

CX6sCX**&. 0 

U«(L*&I*L 

=Ll*fL04TILl> 

uv $3 I*lf 

OC 5j 

AIIw)>C*3 

s: CONTtKU- 

MJ1£3*hUM!3IV-l 

cc u: 

30 i:C 

A3r3<i»j>* }.a 
l:: ccntii.ue 


FCSMCuArS Th£ a kATRIx 

I£U0»MH'C:V»2 
C3 33C :*WfIEFC 

XsX«3X 

CAUL Cfi£f*tC) 

All«1,f3» •♦lt3/3Xo-CUI 
At>3»:«2»«>«.3/cx6«4>a*cm»ci^)»3(n 

I «»l!.3/0Xe-5.3*:(l»-L«a*Ct2»»'2.g*C(3»*0 (4»-v(5» 
AlIoStlM'ZS.OyQXEto.O^C (3)>2.a*C(L)»Ct»l 
AlI-3»I*t) **lS.3/0xe+3.3*C<l»-L.3*C t2l-2. 3»!;i 31 *C (•> *C<51 
Aa«3tl»3) i*6.3/0Xq>4..a*CU> »C(2t »C(3I 
Atl-3,I*3» s*l.3/3Xo*Ctl) 

RMSa-3»«C(Z» 

3:: wCmt;;.u= 

:fCLX«£ .3S. ISa-NU-EC)* CALL N-.?i*.(IPl<lciLlS?» 
w9ir:(S,c) 

33 o3C Isiff^UMSa 
U£riUsMjcSlV*£ 

M-<XrS(6,3) (A<I. J> iJ>ttJ£N3) 
n»ITEIS,ii*» 

LlN£*LlN£,*3 
ft:: CONTINUE 

i»RITEI6»T) 

A-»IT£t6,51 J» ,l»i,KUMe3l 

t.:Ne«LiAE+3 


“OOIFy T>-£ - MATRIX TO ACC3JNT 
f;= icu.-.oARr ccncitiuns 


PAGJ3 fs 
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C»l.k ftCCEFFIBCI 

IFtblKC •6T. CAUL NWPCtlFAfie.wIMEt 

M4ir£(6«UI 

L:NS•Ut^^♦^ 

C!) 3SC Xsl *NUH£3 

left ENO 

A(Xt«)aA(X.«>«A(Xtl>*eCU»2»*A(Xt2»*3S(2*E) 

»:CHT ENC 

m( I*KU(tCXtf »1) «At XiMUHOXV (X* NUJiOIV*^<»)*BCI 3* l> 
l*A(I«t«UF0Xtf»5>*&C(Ufl> 

A(!«Nt;MQXV»2)BAIX»NbN0IV«2l (!• >4U<iOXV»<»)*Bei3f2t 
X*Aa4hU«OXV*S»»8C »Uf2) 

3jA CC.NTXt.UE 


; CISCArO UM-NEE3E0 ELEMENTS 3F ThE MATaXA A. 

JENCaNUOSlV^Z 
O'? 5a; i«i«\uMea 
c; sac jaAiJENO 

m 

&STCCta'«x)<A<x*j) 

« • 

s:: continue 

m 

:F<LXN£ .S£* l5’d-?*MjM£:n CAvL Ni»FUtIFA«»US5l 

-»iT6ie.»T 

00 T3C laltNUHEa 

48.:TH6,3» <A3T3(XtOt«J«l*'iJME'i) 

AFXTE(bfX4l 

hnc«l:l£*i 

7:a CONTINUE 


6ENEFAT6 9 FSCN 4»HS 

IFCLIf.S .3s. r54*NU«eO)J call NWStiUFAUfttLXNsI 
>ialT£l9.9» 

00 5aC I*1»NUM£J 
CC )CC JsltNUMEO 
ia»j>aa«i 
Sta CCNTILOS 

Stt.IiafMSin 

V<^IT£(6,S) 0(X» J) »Jal,NUM£3) 

9CC CCNTINUS 

l:ne*cII'E«>muiie3 


confute EZSCNvAlUSS ANO EISSNVECTORS AN3 ofxnt out E£»Ui.r 
CALL EIC£;ltNOMEQ» 
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C3 TO 33 

ISC') ccNT:nuti 

% - - 

STDP 

■* 

1 PCRTIATt 1H1«2X|*THC £XGErjVAi.ue PROStSI IS PC^HulATEO UITri*t/ */i46H 

t i. KUi*OIV ETA Xf) 

2 PC!«HATt2U9*EF1.0»Al 

I MPHAT(/,13X»*3CWI0ABV CCNOITtOM f4ATU**> 
fC«»ATf3i2Si3.3t/»*2£13*]t 
5 fCBHAT(12£ia.3» 

G f3PiATt/,iaX«*4 ^ATBIX 3£F0B£ BEOiJCTtCN*> 

7 PCB^4TI/,13Xv*«HS*> 

A f’B'tATUi.ax,*^ riATBlA rtlTH 30UN3ABy CCUJITlGriG*) 

J '-3B'*AT(/»2X»*S NATAIX»> 

iW FCSl!*4T(XH 1 

* 

£.*<a 


3U2?.CUTIK£ HVfPG(PAGett*M£1 
:^T£3t~ PACE, LIKE 

Tf*;S SUiBCUTIf.ft iOVA(>CE3 P^I'*TEB TO <•»£<« PAGE 

l:ne«i 

PiGtsPias*! 

APIT£(4tl» PAGE 
1 P3»«iT liHli t2u*» SrtPAGE 1 13) 

S£TtJS}, 

S*.0 




5U9i;CUTIN£ C3EFr<C» 

Th:s nOUINE CALCUtATeS THE CO£FF;>;X£r«TS Of The SCVcfNtAt* 
EOUATION* USErt MUST SUPf<.V STATEHENT3 Cll> THROUGH C(7). 
CO«rtCN/CCEF/ 

»ETA,Xf »X*OX»^tl 

PE4U cir> 

0X2>CX»0X 

:X3»0X2*CX 

OX4aaX2*i:X2 

0x5«0xfc*0x 

c<i)*»n«:/x/2#3/ox5 
:(’)«♦(«. 3-3. a*i«i.i»/t<**2i/cx4 
S(3ta«^U>>a*3-2<»«a*!tLl)/IX**3)/2.C/3X3 
CCs*>«*U.JMLl**2-‘*2.3*‘^tt»72.:) /(X**'*»/3X» 
:<si*»<3.:*»Li*»2-i2.a*5?tu/<x**5»/2.:/jx 
C14»»*«xLl»*3-2. 3*Rti**2>/«X»*5» 

C<7|sSi,1*sT 4/ tETA-l.3>/X-*4 


<5TUM>. 

c'lO 


►RH5INAL PAGE L"' 
il! K)CR QGALi r 


tLeiai^ 



'ju 3 ><)uT:Ne 3 Cceff(!»c> 

ri:s i 3 UTIN£ C£f<t^ 6 ^es THI aOUfOA^V CCf.JlTlEV lafklx, J 3 £^ I'usr 
ioPPLt 5TAT£'«t'4T PUKCTXONS FlHXt«> THrOUOH Pi<*<X.w). 


:c“!iju/co£f/ 

♦ tTAtXr ,X»OK,PU 

sCKi2> 

f'l<»(x,/l}s«-Ca.3-<0/K**2 

F»n<»»os:,a 

><)«♦• e.£/< 

C!x5s3<**"» 

t »,tl xt'DX’PraitTntSUI ♦?.3*P23 <£T^,,-^w1) )/ 

* UX*f 2<*«STi,SU)-2. ;*r?3(ET4,«til I 
S..C2,2»»:,3 

aCUrlMj.S 

3 .:C ?,£J»( 3 <*Pa(*IXF,^tt)-£. 3 *P 23 «XF.^tin / 

«■ nx*? 2 <*tXP,RwlJ» 2 ,i^FZSfXPf?t.ll) 

3CUtX>*< »C(?»l3*l.SI*<2.;*0X2*Fti(iT:,^Ll>-<i.;*Fii(e7A,<win/ 
» (0X*Fl2(€ritRUl-2.J*FtU£TA,r-lM 

* t:<*F12lfcT4,PLl»-2.£*Fil<e:TA,^Ui)» 
a; H,2l*'j::«*Fl2l£T».tnUl* ♦2.3»Fit tir*,iirt.U t/ 

► (3X»Pl2«ETA,ffLl>-2..'*?lHiTA,R.|.i>l 

sC U«i> sOX*Fi2<sF*=U)- 2*3»^U{ <Fi -til > / 

♦ i3.<*F;2(*?,ui»»2.o»FitUF,aLi>r 

*■ (2.3*FtllXF,RLtM'CX*Fl2(xF,r!.in 

•-Mv 13 , 2 )-l. 3 »*.» 2 . 3 *’ 3 X*Fi 2 {x?, 5 !i.l)-‘:< 3 *=’tVUf ,RLi)>/ 

*■ (.:**Fl 2 txF,XU»» 2 , 3 *FmxF,^tin 

'40 
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SU3itCUrtNS EXG£KINUKeai 

THIS 3U3R0UTtT.E CALLS THE IHSl EOXTION 5 ACUTIr^E EICIF T3 IC.iPUTE 
The SlCENVALLiTS AnO kXCEHVCSTO^S OF HAT^XX A A*tO liATAXX 3. 

CCHMCK/EXCCMW/ ' 

*A(33*3$)«B{3S*39i 

^EAL a, 3 

COMPLEX 

*ALfA(3S)tLA>l30AI351 •Z(35»35) 

seal 

»aETA(33l«3K(2E0a) 

AaiT5(4«l) 

K^XT:t6.2l NUMEJ 

CALL ElCZP(A«33t3«3StNOMS2*2» AL*At3ETA,Z,A5«/«Xt EER) 

.lt* i.o) FSFRi-iaHwsLL 

:F(ri<ut .ar. i.s .anc* ukcu .lT« ljz.h p«frh«:£hsatsf<*ct.v 

IFr<i*<i» .5T, tOw.61 PRrSf'atOHOflOSLT 
h»ITE(%«6l PRF5.H 

wRITs(o,l3» lER . 

OC 23C Isl.SUWEJ 
IZv LA<<sOA(:t«ALFAi:i/liETA(:i 
hs;T£lc«7l 

rtRITcta.dl CLAH3CA»I>,I«l,NUM5a» 
i*RITEiat^» 

SC 3 3: intNUt^EC 

»RIT£(e*A) IZII«J)*JalfNU.lEC) 

I3C i<R:TE(6,3) 

1 F:RHAniHl»39H Tiid SISENVAlOE *R09fc£M *J SCwrfSC '♦ITHf/t/.l-H 
♦■MJHEO 

2 rCPHATJUil 

I fsrmahih > 

e FSBHAT UmC»34M NOTE— THE ALGORITHM h«3 xCB<«£ 3 |Ai:i 
r FCe -f^lTt j.NC,;»3rt CICENVALUES LAk33L(1) L»H-33A UiUHECt AR£ ) 

? F09HAT(<»(iHit£l2«6«AHt,£12.6«liMI 

9 FCRHATtlHCtJSH SISENVECTSRS AP£ SUCCESI/E COLJInS OF) 

13 FCSHATt IHOtlOXtoMlcr. s ^13) 

RtTUSN 

END 


POOR quality 


TH£ € XGf \u rOK#*yC-T£D 41 tH 

t €T 4 xF 

^ U **>2^Q X.QumS 


£ SiT-lpt 


«^L ^-^UCIXCK 



C* 

■ t ?■*» t ♦ B a 

-• 4bi*£t 0 9 

•&7DC*99 


» 7 juL*Cd 

i« Q« a* 

G* 


c. 

d. 


£ 4 

..?J5E»09 

t. 

-.^Z«.£»C9 

«:»90E*a9 

■** JJ2£*C9 

»7^nk43« 1]. 

4* 

a. 

0* 

a« 

t * 

4<.?3t»B7 

-pb9uC«fl d 

»* • 

•2aCE»99 


■ 6 J.J£*S^ 

-•J2o£*C9 .7iu£«}4 y« 

3* 

a. 

a* 

a* 

C« 
> * 

V « 

C. 

*762t^H7 

• • 

-•d2>£*0!l 



«bll£4Ca «^.32J£»39 *641£^B« 

3« 

a. 

c. 

I* 

C * 

z* 

c . 

h)« 
■# • 

■ SC3k*aci 

^•^57c*B4 

4iCl,£*D9 

^4u79t4ii .bioEta^ *«3tSc»a9 

•«97c»4t 

y# 

a, 

s* 

z* 

3. 

5 * 


tJ • 

• I23£4^34 

*.tj7d*59 

•3otC«ii9 *»ob4i£*39 

49 

• ti3bE*34 

a* 

«6ie€«:rt 

u 

-« • 


« 
^ 4 

ft* 

3. 


-•U7E*S9 .4»9£*39 

•E21fct44 

*4»«yC*49 

-• 45Jt *3i 

w « 

40- 4^ *3j 

t> * 

C 4 

c. 

d4 

a. 

V* 

•197^*90 *#ii:**E»39 •btie*:9 


•*7J£*49 

B. 

4* 

-■•5 ja-L *C'f 

4£*C0 

V • 

w • 

M « 

a* 

y • 

*173E*a4 ^■132£*G9 

•4tlSC«39 

-•b4JE*a9 

a. 

B* 

• t 

J * 

^ 4 2-j5t4BTl 

34 

• 5^77Ftae 

a* 

a. 

G. 

9* 3. 


•b2«£*49 


MS 

••*^^» c*i *3 1 

<ic+:i 


*.3X4£*B2 

*#i7sc*a? 

-•95SE*C1 

*,5$7£»4x -.iji»e*ai «'.269C»at 




^CUNJ*.-* ? * 3 r*aiTI 0 h 

9 » 

• l03Ei-9l *b3mw*^*51 













« 7*f U&4 a 8 

ft* 

^* 

ft* 

tf* 


a. 







. 7^bir^ C A 

ft* 

s« 

ft* 

3* 

ft* 

1 

i 

^w«S<E«3>t 

* 


• 6-i3EHJ5 

-*32iVt*ft'si 

*71ftt*0fl 

0* 

ft* 

ft. 

ft* 


f 




-•6ft;c»cs 

*bltc^ft9 

**32ftE»69 

*<MSlt4 08 


ft. 

0* 



3* 

-i:3L#ft« 



*.62&E*3*3 

,.61bE*C^ 


*o&7E*ftA 

ft. 

ft* 



3. 

3« 





•bl 



a* 



<>* 

a. 


.i4>iE*aa 

-*117E*34 


**€ii ftc>^ C9 

*b£*E «-ftS 

*.i3«£fft9 




a. 

4i« 

0- 

u ■ 

*i&7£«aA 

-.*^^£♦39 


-*o72tt*^f 

*b£2£*ft9 




0* 

y« 

» ♦ 

ft. 

ft* 

•i73C«a8 



<^*oi»ft€:«a9 

*5«9E«d9 

'sj 



a« 

3« 

a • 

ft* 

ft* 

ft* 

•143£*B4 

**l3<*Et09 


**6ft«e*89 

O 


u f*Ar-^i< 












-*73Zf ^13 

4 • 

- • 

; « 

ft* 

ft* 

ft* 

ft. 

ft. 

a* 



3, 


c « 

ft » 

ft • 


ft* 

ft* 

ft* 

w* 


j 

&. 

#• 

-, :31Li53 

V ■ 

ft* 

ft* 

■j* 

I.* 

ft. 

ft* 



c* 

y« 

a* 

i£*lJ 2 

ft* 

ft* 

ft. 


ft* 

ft* 



c* 

J* 

y 4 

* « 


c. 

41* 

J* 

ft* 

ft* 


{ 1 

d* 

J» • 

G# 

* 

V . 

9* 

-*170£^ft2 

ft. 

ft. 

ft* 

ft* 


i J 

d* 

s* 

4 • . 

t 

a* 

ft* 

£*Ci 

*. 

J. 

ft* 


^ i 

ii. 

c. 

3. 

y * 

ft* 

ft* 

a* 

-*5ft7t*ftt 

ft* 

ft* 


i 

a. 


3* 

ft m 

0* 

ft* 

ft* 

ft* 


ft* 



ft# 

y* 


• 

ft* 

6* 

ft. 

0* 

A* 

«*269£»itl 


J 

1 

i 











tMC P«>eu« rs SCLVtO UlTh 

•^U*-fcC 

n 

liCrf— IK. itt;c.;:THrt hA;; wO'<SJ -<ECL 
lEK * C 

r b*. jJitfijMlii ii<£ 


1 

II p^CSr^nL^-M jp«t p 
M J7?f *37,3. 

1 ( p<«9<p7: 

It p242*«ei£*u6T0* 

I f » J1 j^^buC *-unpb • 

I I 

IE .6 J^0«kQ£*04|C. 
1 ( ptl7<*53(^«^37,Q. 


t !(!(,?-* r *; ii<£ 

f 5^ J £ *V“* m'* f ^ « 

1 pt S€?rw? *-i3t 4# 

jU3C£>icr ^s Cl 

> l b7»ti£~3£vSp 

1 1*. 27e>^<«^£*^Q0«mp 

II pl25li7£-Ca,3* 
It plCC035C«i:i*0« 
It 

M-.3 34*2J3£*03«d, 
It-P f J*p72lpcb33«C* 


1 **vlt(* 

1 « t SC » a « ' * M I 1 j * 

1 p A f S G ivr ♦ C3 t ’ • 

M p 

kt pp4 lcil,t>*3C*v p 
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APPENDIX II! 


Free Surface Boundary Conditions 


Free-surface boundary conditions on temperature are Identical to 
rigid surface boundary conditions for temperature. Velocity boundary 
conditions for rigid surfaces are (Durney, 1968a) 

P|_ = Pl" = 0 at r = n t 1 . 

From equation 10(c) we find the equivalent conditions 


and 


Dj^ 6^ » 0 at r = n» 1 


[r^ D, e, ] * 0 at r = n » 1 
dr^ L L 


By Inserting the definition of the operator these conditions can be 
written 


d'* 


+ 


8 ^ 
r ^ 


dr“ 


+ 


(14 - L,) d^e L 


dr" 


4-2L, d0, 

+ — ^ ~ *0 at r * Tii 1 

r"^ dr 


and 


d^e. 


. 2 


d0| L-i 

— ~ ~ “t = 0 at r = n» 1. 


dr dr r' 

The mode of the free surface solution for wave number L 
shown in Figure 11 . 


= 3 is 





free surface 


Figure Ub: Free Surface Solution: Temperature Profile. Dotted Tirie is 
temperature profile of pure conduction, solid line shows sum of eigen- 
vector 0j^, L = 3, and conduction profile, has been arbitrarily 
normalized so its length is unity. 
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THERMO* HYDRO-DYNAMIC CHARACTERISTICS OF A • 

( 

KRO-ORAVlTYj SPHERICAL MODEL OF THE TROPOSPHERE | 


S. SrlvatBuign 

•- j , RcMBrch Assocl«t« „ ■*>. f .authors 

D«pt. of AtM*ph«rtc Sclmco 
) CeloraAo St«to Unlvorticy 
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ABSTRACT 

j 

A aodel that eaptoita tho radial irtartia forcaa of 
a rotatiitf fluid containad in a spherical annulus 
is deserifaad. The nodal would ba flown in a satel- 
lite and exparinents would ba perfomad in very low 
or taro gravity. In such a nodal it would not ba 
nacassary to artificially tinulata a radial gravity 
fiald. thus snail anounts of alaetrlcal energy 
would ba sufficient to perform exparinents. Since 
the only forces involved are thermo- h/dro-dynanic 
ones, alectronagnetic aquations need not ba consid- 
erad. 

INTRODUCTION 

A variety of experlaants have bean perforsiad under 
usual laboratory conditions to simulate the large- 
scale dynamics of the earth's atmosphere^'*''. Of 
these the most successful so far are the cylindri- 
cal annulus experiments in which a liquid is con=. 
fined between two concentric circular cylindars^*^. 
A radial thermal gradient — resembling the meridi- 
onal temperature variation in the atmosphere --is 
imposed. Cameras are attached to record notions at 
different levels in the fluid^^i. The entire 
system, including the cemeras, is rotated in the 
same sense as the earth. Thus all observations arm 
made relative to the solid-rotation rate. fl. of 
the system, and they resemble observations of the 
earth made by geosynchronous satellites. 


fora these experiments do not possess an important 
characteristic of the earth-atmosphere system. 

A spherical annulus model (see Fig. 1). which has 
meridional variations of temperature and the local 


9 

o 



Fig, 1. The ei^riaental setup for the taro- 
gravity, spherical atmospheric model. 


The cylindrical annulus experiments have some of 
the basic properties of the earth- atmosphere system, 
vit., a meridional thermal contrast called baro- 
clinicity, and a non-Newtonian frame which is due 
.to the solid rotation. Under laboratory conditions 
both these properties can be varied. Such varia- 
tions have led to important results concerning the 
breakdown of toroidal (or Hadley) cells, the estab- 
lishment of wave regimes — which are asymmetric 
with respect to the axis of rotation — and the re- 
establishment of the toroids^*'. 

It may be noted that the tropics are dominated by 
toroids (with a *h«eah** wave regime) and the extra- 
tropics by wave regimss"f. Thus axially syamietrle 
and asymmetric regimes coexist on the earth. This 
coexistence is due to the meridional variation of 
.the local normal component of D. In the cylindri- 
cal annulus models the angular velocity D is eueiv* 
.where., normal to the base of the cylinders. Theret_ 
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normal componont of cannot bo used under normal 
earth gravity to simulate tho bohavior of the 
atmosphttee Such a model may be used successfully 
under very lov or sero gravity conditlonSe It has 
been suggested that a zero-gravity, spherical model 
must hava a radial gravitational field* simulated 
by the imposition of an electrical force on the 
fluids The consus^tion of electricity for this 
pilrpostp howeverp Is a significant fraction of chat 
available in a satellitt^**^* Hence an alternative 
should be founds We suggest here that the inertia 
forces of the spherical annulus model may be used 
to simulate the effects of a radial gravitational 
field. 


.ilBCf.r.i"''- " ■ nn W>WK NOT TTX- 



• ?>- 


tf.r P.’ 'J . 11 P**- 


'FORCES ACTING ON THE FLUID IN THE HODEL 

It is wcll'knoim fro* th* governing equetionf.of 
■eteorology [see, e.g., Haltlncri”, Holton^) 
that the forces to be considered in neteoroiogy are 
t. the inertia forces due to the rotation of the 
coordinate systen, vis. the centrifugal and Corio* 
Its forces, 11. the prcssure«gradlent forces, ill. 
the buoyancy force, Iv. the viscous forces and v. 
the gravitational force. Other kinds of forces are 
of negligible inportance in awteorology. 

No note that (static) atnospherie pressure is due 
to gravity. Thus pressure-gradient forces are 
Indirectly due to gravity. 

In Fig. 2 we present all the possible coaponents of 
the inertia forces in a spherical, coordinate systeii. 

a 

NORTH POLE 


r fg / 


EQUATORIAL PLANE 


SOUTH POLE 


Fig. 2. Ihe coaponents of the inertia forces in the 
spherical nodal. The synbols. are described 
iy Eqs. 1 to S. 


Ihese are represented nathenaticaUy as follows: 

• rfl^ eos^g, 

(radial cantrlfugal force),. (1) 

• -rO* cos ♦ sin d, 

(north-south centrifugal force), (2) 

t. • 2Ihi cos d, 

(radial Coriolis force), (3) 

j f. ■ -20u sin d, 

) (north- south Coriolis force), (4) 


and f, • 20v sin d ' 20w cos d» 

• # a 

*3,1 *3,2 I 

(east-west Coriolis force). ($) 

In the above equations r is the radial distance 
fron the eoasson center of the two spheres In Fig. >, 
fl the angular velocity of the systen, d latitude "'* 
angle, and u, v, w the eastward, northward and 
upward components of motion, respectively. The 
coaponents of motion are positive if in the sense 
mentioned above. 

4 

The distributions of these forces clearly obey the 
following limiting conditions: 

»j, (i»j, fj, fj ^2 0, as d ly (*) 

« 2 * ^ 2 * ^ 3,1 ♦ 0 * * 0 

Combining Cqs. 1 and 3 we obtain the radial inertia 
force equation: 

2 2 

gy • rll cos d *■ 2(hi cos d (d) 

t 

In Eq. I the left-hand side term la analogous to 
gravity and hence is denoted by g; the subscript H 
stands for ''model'*. Let us determine the relative 
magnitudes of the two right-hand side terms of Eq.8. 
The ratio of these two terms is: 


2Du cos d 
rO^ coa^d 


2u 

Or cos 4 ' 


The right-hand side (rhs) of Eq. 9 represents a 
Hossby number, denoted by Ro. 

It is well known in meteorology that the condition 
for the prevalence of quasi-geostrophic equilibrium 
is that Ro « 1 (see, e.g. , Holton ^°'j. Since 
tropospheric notions are quasi-geostrophic, and 
since we want to reproduce and study such conditions 
in our model, we may assume that for our model 
arperiments also Ro « 1, (The validity of this 
assumption can be experimentally established by mak- 
ing small. See Eq. 16 below.] Using this value 

in Eq. 9 we see that: 

30u cos d « ces^d. (10) 

Thus the second rhs term of Eq. 8 is negligible 
against the first rhs term. Therefore Eq. 8 
reduces to: 

• rfl* cos*d. (11) 

Equation 11 may be interpreted as follows: In a 

spherical, inviscld. homogeneous, rotatine f^utid mass 
not scyd OT by gr.ivttationai forces, the rndlal 
centrifucal force .icts as a spherically asymmetric 
gravity-like force. 


» t 
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CHARACTERISTIC PRESSURE DISTRIBUTION IN THE 
• WDEL 

1h« hydrostatic equation ia: 

dp m • g p dr» (H) 

whoro 6p/dr ia the variation of prciture with the 
(verticai) coordinate r, g the acceleration due to 
gravity, and p the density of the fluid. Ut ua 
assuM, for the aonent, that p is constant. Substi> 
tuting Eq. 11 in Eq. 12, and integrating froa r, 
the radial position at which the pressure is to be 
deterained. to r,, the radius of the inner sphere, 
we obtain: * 

P(t» ■ I 0* (r* - rj[) cos^q. (IS) 

Froa Eq« 13 it aty bo ooon that the pressure duo to 
the radial centrifugal force is a oaxiiiuBi at the 
outer sphere radius (r • r,} at the equator 
C4 * Op cos ^ * l)e The pressure dininishes to 
tero at both poles; it is also tero at the inner 
sphere (see Fig. 3)« 



Fife 3* The distribution of static pressure due to 
the radial centrifugal force in a spheri- 
^ilp tero^^gravity model of the atmosphere. 

Units: dynes cm*^. 


DISTRIBUTION OF ZONAL NOTION DUE TO INERTIA 
FORCES 

From jcq* 13 the equation for the north^south pres* 
sure*gradlent force is obtained as 

' p r^ • 4 

Conbirlng thia fore# with the other north~south 
inertia forcta (Eqs. 2 and 4) we obtain 

^ ~ (r* - rj) sin * cos d - 

* cot 4 sin d - 20u sin b* (1$) 


/ X ^ I fKt :'■*» 

If the meridional acceleration dv/dt f Op Eq. IS 

u ■ .joXees#. I (W) 

the following points concerning Eq. 16 ere worth 
noting: 

t. At the equator (6 • 0) every one of the rht 
tenu of Eq. 15 it aero. Ilnis u is indeteninAte 
at the equator. Therefore, Eq. 16 is invalid in 
the vicinity of the equator, and tht naxlnu* aagni* 
tudaa of tht zonal notion due to inertia forces nay 
ba e:qpacted in the niddle latitudes. 

11, la the above coosiderations, tha effects of 
viscosity have not bean taken into account. For 
any viscous fluid a no-slip condition has to be 
applied at tha spherical walls. Hence thd fluid in 
contact with the welts will hsve only the solid- 
rotation velocity. 

ill. Since the enstcrlies represented by Eq. 16 
are dissipated by viscosity at the walls, angular 
wmentisi wlanee in the aodel requires that wester- 
Itea be dissipated elsewhere. This aight be 
expected to happen in the tropics, where viscous 
forces eight help.estebliah toroidal calls in the 
Banner of Pearson^*'. ’ 

iv. Since easterlies are established in the 
■idletitudes by the inertie forces in the spherical 
annulus nodal. It is essential that tha nodal be 
heated In the equatorial regions, and cooled near 
the poles to produce westerlies In the aidlati- 
tudss, as in the earth's ntaosphere. We also not# 
froa Eq. 16 that tha oastarlias established by 
inertia forcta zan bt reduced by decreasing the 
value of rj^. If Tj > 0 (i.e,, for a spherical 

fluid aasa), tha zonal notion due to inertia forces 
is identically zero. 

GENERATION OP THERKAL CONVECTION 

It was shown earlier that the radial centrifugal 
force la nuch larger than tha radial Cori-'lis force 
under qunsl-geostrophic conditions, and that the 
radial centrifugal force is a function of radial 
dlstanea r and latitude 6> When the radial centrif- 
ugal force is tha doainant radial forrr the wara* 
ing of tha equatorial regions near rj, inner 

radius, will not Itsd to theraal convection, since 
only denser fluid is drawn to the outer radius r,. 
^erefore, to siaulye cellular conyecttc r . the * 
fluid aust be heated at Tj rather than at ' r j . We 

note also the following: 1. Even if the entire 
outer sphere is heated, convection will nut be 
generated nt the polea, since the radlnl centrifu- 
gnl force ia zero at thasa points. 2. The heating 
of tha outar aphere leads to a radial distribution 
of taapornturo in tha aodel corresponding to the 
varticnl (redial) distribution of potential cenpera- 
tura in tha ataosphera. This is the proper siai- 
larity between the inconpretslble fluid in the 
aodel and the conqirestibla ataosphera. 
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’CONCLUSION 

It c«n bt rtHilUx shoim that tK« aeridlonal pr*s- 
sur# sradient due to th« redial cantrifugal force 
(Eq. U) leads to the establishment of westerlies 
which increase in magnitude radial ly» under geo* 
strophic c^iiditions. However* the meridional 
centrifugal force counteracts this, and establishes 
easterlies in the middle latitudes (see Eq* 16)* 
Therefore the establishment of westerlies in the 
midlatitudes of the spherical annulus model demands 
the imposition of a meridional temperature gradi* 
ent* It has been shown above that heating the 
outer sphere* rather than the inner sphere* leads 
to thermal convection under geos trophic conditions* 
If* however* zonal velocity u becomes comparable to 
fir* so that fio 1 (see Eq* 9)* the radial Coriolis 
force also becomes important in establishing 
thermal convection. Under such conditions the 
equations of motion are quite non-linear* and the 
resulting circulations have to be studied through 
actual experiments or through numerical models* We 
have begun a numerical model study at this time. 
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ABSTRACT 

Very good agreement !a shown to exist between the meridional distributions of the sonel wavenumber n of 
rapidly amplifying baroclinic waves on a sphere and of an average wavenumber It of “grid'Sede” atmo> 
spheric eddies. As a consequence, the xonal wavelength both baroclinic and atmospheric eddies remains 
virtually constant, i.e., within a factor of 2, over the extrat^ics. The values of m at different latitudca 
have bm obtained by using linearised baroclinic theory on different meridional profiles of the unp^turbed 
sonal wind (MPUZW). Since they agree with V, atmospheric eddies are, in relation to linear baroclinic waves, 
independent of MPUZW. In tbit tense. V is controlled locally rather than globally. 

The mutual dependence of the upward and poleward transports of (Knrible) heat in baroclinic wave thcMy 
is cottecUy formulated— as compared to a direct analysis of the first law of thermodynamics. 


1. Introduction 

Mount and Stone (1976, hereafter referred to as 
M ft S) and Stone (1974) have derived a number of 
interesting properties of linear baroclinic waves. They 
have compar^ their results with available meteorolcig* 
ical statistics such as those of Oort and Rasmusson 
(1971) and found quite good agreement. In the case of 
two specific results, additional methods can be used to 
infer the properties of atmospheric eddies, which may 
then be compared with the M ft S results. The methocU 
we have alluded to will be discussed in this paper. 
The results under consideration are the zonal scale of 
eddies and the interdependence of the poleward and 
upward transports of hrat 
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Fro. I. The latitude dependence of xonal wavenumber for linear 
batocUnk waves oa a sphere; dots, near-neutral itability (aym- 
metric itieamfuncsiao); aided dots, far from neutral stabiUty 
(f)'tnmetiic lUeanAnic^); dieted crow, ftr from neutral stab- 
ility (streamfunctian antiaymmctric with cexpect to the equate:). 
Abacusax of points refer to the totitudes at which the geopotential 
dgcnfunction, at the top level of the threc-lcvd model of Moura 
and Stone (1976), reachet the maximum magnitude. Each point 
conexponds to a different metidional profile of the unperturbed 
aonsi .\ind. 


2. The zonal scale of eddies 

For the study of this property we conrider the 
spherical, two-layer (three-level), linear quasi-geo- 
strophic model of M ft S. M ft S use several meridional 
profiles of the unperturbed zonal mnd (MPUZW). 
Their results concerning the zonal scale of linear 
baroclinic waves may be summarized as follows: 

1) The nearer the peak in perturbation geopotential 
is to the pole, the smaller is die zonal wavenumber (is) 
of the most unstable mode (see Fig. 1). 

2) As a consequence of 1), the zonal wavelength near 
the peak in perturbation geopotential is about 5000 km 
on a sphere of the size of the earth. (This wavelength is 
virtually independent of latitude in u broad mid- 
latitude belt) 

3) The zonal wavelength is pn^rtiona! to the 
radius of deformation at the latitude where the max- 
imum of perturbation geopotential occurs. 

Hie author (Srivatsangam, 1976a, b) has described a 
parametric method which yidds an average zonal 
wavenumber U for atmospheric eddies throu^ the 
geostropbic meridional wind equation. Herei, the 
overbar detiotes a weighted toot-mean-square averaging 
(for details, see Srivatsangam, op. cU.) Now, 

»-(C***K«**])*(rV«cosd), (1) 

where square brackets denote zonal averages and 
asterisks departures therefrom. Hence 
represent the zonal variances of the geostrophk 
meridional wind at, and the bright of, an isobaric 
surface, respectively. Also, / is the Coriolis parameter, 
a the mean radius of the earth, p latitude, and g the 
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Tabu 1. Height and latitude variatiou of the monthly mean value* If (ID and the temtMtal coefRcicnt of variation Cfftf 

of the average aosal wavi aumber 00- 


PrettUft 

(mb) 

20 

25 

30 

35 

40 

45 

Utiuidt CN) 
SO 55 

60 

65 

70 

75 

80 

85 

100 

4,02 

3,51 

3.77 

3,87 

«. 

4.10 

If (R) for October 1968 
3.92 3.29 2.70 

2.25 

1,86 

l,S9 

1,39 

1.27 

1.11 

200 

423 

4,26 

4,45 

4.71 

4.83 

4,57 

4 .03 , 

3.47 

2,93 

2,41 

2j05 

1,79 

1,49 

U5 

300 

5.13 

4.8S 

4.95 

5.16 

5.10 

4.76 

4,27 

3.77 

3.25 

2.72 

2,36 

1,98 

1.62 

1.32 

500 

5.61 

5,31 

5,26 

5.44 

5.03 

4,64 

4,20 

3.62 

3.11 

2.6S 

2J7 

1.89 

1.61 

1.27 

700 

5,05 

5,17 

5.24 

5.14 

4,84 

4.53 

4.03 

3,49 

3,07 

2,63 

2.21 

1,88 

1.61 

1,24 

1000 

4,95 

4,98 

5,26 

5.30 

5.16 

4,68 

4,08 

3.63 

3.35 

3,07 

2.70 

2,09 

1.55 

1,14 
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acceleration due to earth’s gravity. From (1) it is seen 
at once that H represents "grid'Scale” atinospheric 
eddies. 

The daily values of H were computed at the 100, 200, 
300, 500, 700 and 1000 mb levels in the region 20°N to 
SS^N using the National MeteoroIogical'Center (NMC) 
data for October through December, 1968, and February 
through April, 1969. The hour of observation of the 
data used here was 1200 GMT, and a 5” latitude by 
5” longitude grid was used. 

The data on monthly mean values of lit, denoted by 
If (H), are presented in Table 1. These are for October 
1968 and February 1969 only, but are typical of all 
the six months for which computations were made. 

From Table 1 and Fig. 1 it is readily seen that result 
1), and hence result 2), of M & S are in good agreement 
with the ensemble characteristics >of atmospheric 
eddies as represented by It. 

It may be repeated here that each data point in 
Fig. 1 corre^nds to a diftrttti MPUZW. However, the 
data of Fig. 1 and Table 1 show excellent correspond- 
ence. From these, the following conclusions may be 
made: 

In order to derive the typical wavenumber of 

PiVlrK ^ 


atmospheric eddies in any particular latitude through 
linear baroclinic wave theory, a MPUZW which would 
yield a maximum value of perturbation ge}potvntial 
height at that latitude must be used. 

(Conversely, the zonal wavenumber of atnosphcric 
eddies at ea^ latitude is such as if the maximum 
ma^tude of perturbation gcopotential height occurs 
at that latitude. 

Also, from the last statement and result 3} of M & S, 
it follows that the zonal wavelength of atmospheric 
eddies is proportional to the radius of deformation at 
each latitude. [This does not contradict result 2) of 
M S. Since the radius of deformation is proportional 
to 1//, which varies by a factor of 2 from 3U°N to 
90*’N, the typical wavel«mgth of linear baroclinic and 
atmospheric eddies varies only by a factor of 2 in the 
extratropics. On the other hand, the wavenumber, os 
seen from Table 1, varies by a factor of 5 from the 
subtropics to the subpolar region. By comparison, 
therefore, the wavelength rather than the wavenumber 
is constant across latitude circles.] 

Now* wC turn to thequestion of the temporal variiibil- 
ity of It. Table 1 contains data On the temporal co- 
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Tabu 2. Height and latitude variationa of the temporal coeOicienta of variation C(£A*'J) and C(£«**3) Ihe tonal 
variances of itobaric surface height and geosliophic meridional wind. 
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efficient of variation of 3 given by 

o(S) 

C0l)-_X100, (2) 

Jf(K) 

where (t( 10 is the monthly standard deviation of If. 
Using the data described above, C(S) was calculated 
for all latitudes, pressure levels and months. 

The C(8) data in Table 1 indicate rather small 
day-to-day variations in m. A comparison of the 
coefficients of variation of C***3 (see Table 2) 
and If shows that C(il} is only a^ut one-half as large 
os C(Qs**3) and C(Cv**]). From (1), then, it follows 
that variations in Cs**3 are compensate by like 
variations in This is confirme by Fig. 2, in 
which the daily values of C»**3» C*'”] ^ 6S®N, . 

200 mb and October 1968 are presente. 

The temporal quasi-constancy of H thus achieve by 
the mutual compensation of [s**] and [v**3 indicates - 
that the results of M & S are valid not only on a monthly 
mean basis, but also on a daily basis, i.e., for transient 
atmosphaic stater. 



Fic. 2. Daily values of the aonal variances of the geostrophic 
meridional wind Ct*'] and the height of the 200 mb surface 
and the weighted iooi-mean>tquare wavenumber R. Values are for 
200 mb, 65*N and October 1968. 


Finally, the values for October 1968 and 

February 1969 in Table 1 reveal rather small inter- 
monthly changes. This is also true of the other months 
for which computations were made. Thus it appears 
that in each zone ?! remains nearly constant at least 
through the cold half of the year. This may be compared 
with the phenomenon of vacillation in the cylindrical 
annulus experiments (see, e.g., Lorenz, 1963, and 
Pfeffer and Chiang, 1967). Vacillation is a process in 
which, under constant conditions of rotation and 
impost radial thermal gradient, the wavenumber 
remains a constant whereas the wave amplitude or 
radial axis tilt varies cyclically in time. Since « remains 
nearly constant in the atmosphere, a similar ensemble 
average wavenumber may perhaps also remain con- 
stant in unsteady Rossby regimes. These arguments 
are, however, quite conjectural, and need verification 
through d.ita analyses. 

3. The interdependence of vertical and meridional 
eddy heat fluxes 

Stone (1974) has derived the following equation foi 
the upward flux of the sum of sensible heat plus 
potential energy Qvhich may be closely appro.ximated 
by the upward flux of sensible heat (see Oort and 
Rasmusson, 1971, pp. 50-51)3: 

C6*«,*3=-2»c<lir*,-rM*-}-2C«3C»*»*3- (3) 

In (3), 9* is perturbation potential temperature, u, e, 
u> the zonal, meridional and vertical components of 
velocity, respectively, «•**, #*i the amplihute of perturba- 
tion pressure at the top and l^ttom levels, respectively, 
and c< the imaginary part of the complex phase speed. 
£q. (3) applies to baroclinic waves on a plane. 
Since CjeO near neutral stability, D9*w*3 >* directly 
proportional to C0*v*3. M & S have extended this result 
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of Stone (1974) to barodinic states far from neutral 
stability, and on a sphere, and concluded that even for 
these cases is proportional to Therefore, 

the first right-hand side term of (3) is negligible as 
compared to the second. 

Kuo (1956), examining barodinic instability in a 
cylindrical coordinate system, has derived equations 
similar to (3) above CKuo, op. cit . ; Eqs. (72) and (73)]. 
The numerical values givep by Kuo for the upward flux 
of sensible heat are in order-of-niagnitude agreement 
with meteorological observations. 

Srivatsangam (1976a) has derived the following 
equation for the upward flux of sensible heat, from 
the first law of thermodynamics: 

= -l»»{ra-*|dCr«]/2a/-i-Ct-](3Cr**]/2fl3*) 
+Cw](dCr*]/2a{)-i-C»T*]aCr]/(id^| 
-}-4-*.lCr*a«*/dO+C«][r*d9Vff cos^x] 
+Ct-Xr*39*/fla^]+C«']tr*ag*/3{] 

+E«*T*](dC93/a^*)l 1 -5* (H-y)- (4) 

- In (4), Cp is specific heat of air at constant pressure^ p air 
density, T temperature, the dry adiabatic lapse rate, 
being the latent heat of vaporization of 
water, q the specific humidity of air, 0, X, ( and / 
latitude, longitude, altitude and >.ime, respectively, and 

r - rr»a[r]/d{+ 

Eq. (4) involves several assumptions including the 
omission of diabatic effects other than the release of 
latent heat of vaporization of water and the negligibility 
of triple correlation terms such as [X*v*dT*/aB^2* 

We see that only two terms on the right*hand side 
of (4) involve the factor An order-of-magnitudc 

estimate of all right-hand side terms of (4), using the 
Oort and Rasmusson (1971) data as reference, shows 
that the term involving is at least 

one order of magnitude larger than the other terms. 


Thus the formulation of the interdependence between 
upward and poleward fluxes of sensible heat in barcK linic 
wave theory is in agreement with a direct analysis of 
the first law of thermodynamics. However, very many 
other effects are neglected in deriving an equation like 
(3), although these effects are small Thus, order of- 
magnitude agreement between the predictions of 
linear barodinic wave theory and meteorological data 
may be expected even in the case of upw ard eddy heat 
flux. This is indeed proved by the restdls of Kuo, 
mentioned above* 
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